Working out the bugs in my new object that implements assembly bias in a continuous way.


In [1]:
from pearce.mocks.kittens import cat_dict
import numpy as np
from scipy.stats import binned_statistic, linregress

In [2]:
from matplotlib import pyplot as plt
%matplotlib inline
import seaborn as sns
sns.set()

In [3]:
sns.set_context("notebook", font_scale=1.5, rc={"lines.linewidth": 2.5})

In [4]:
cosmo_params = {'simname':'chinchilla', 'Lbox':400.0, 'scale_factors':[0.658, 1.0]}

In [5]:
cat = cat_dict[cosmo_params['simname']](**cosmo_params)#construct the specified catalog!

In [6]:
cat.load(1.0, HOD='redMagic')#, hod_kwargs = {'sec_haloprop_key':'halo_log_nfw_conc'})#, hod_kwargs={'split': 0.5})

In [7]:
cat.model.param_dict


Out[7]:
{'alpha': 1.02,
 'f_c': 0.19,
 'logM0': 12.2,
 'logM1': 13.7,
 'logMmin': 12.1,
 'sigma_logM': 0.46}

In [8]:
fiducial_point = {'logM0': 12.20, 'logM1': 13.7, 'alpha': 1.02,
                      'logMmin': 12.1, 'f_c': 0.19, 'sigma_logM': 0.46}

In [9]:
rp_bins = np.logspace(-1,1.25,15)
pi_max = 40

In [10]:
hod_params = dict(fiducial_point)

In [11]:
n_pops = 2000
wp_vals = np.zeros((n_pops, rp_bins.shape[0]-1))

for i in xrange(n_pops):
    if i%10==0:
        print i
    cat.populate(hod_params)
    wp_vals[i,:] = cat.calc_wp(rp_bins, pi_max, n_cores = 1)


0
10
20
30
40
50
60
70
80
90
100
110
120
130
140
150
160
170
180
190
200
210
220
230
240
250
260
270
280
290
300
310
320
330
340
350
360
370
380
390
400
410
420
430
440
450
460
470
480
490
500
510
520
530
540
550
560
570
580
590
600
610
620
630
640
650
660
670
680
690
700
710
720
730
740
750
760
770
780
790
800
810
820
830
840
850
860
870
880
890
900
910
920
930
940
950
960
970
980
990
1000
1010
1020
1030
1040
1050
1060
1070
1080
1090
1100
1110
1120
1130
1140
1150
1160
1170
1180
1190
1200
1210
1220
1230
1240
1250
1260
1270
1280
1290
1300
1310
1320
1330
1340
1350
1360
1370
1380
1390
1400
1410
1420
1430
1440
1450
1460
1470
1480
1490
1500
1510
1520
1530
1540
1550
1560
1570
1580
1590
1600
1610
1620
1630
1640
1650
1660
1670
1680
1690
1700
1710
1720
1730
1740
1750
1760
1770
1780
1790
1800
1810
1820
1830
1840
1850
1860
1870
1880
1890
1900
1910
1920
1930
1940
1950
1960
1970
1980
1990

In [12]:
sns.set_style("darkgrid", {"axes.facecolor": "0.85"})

In [13]:
rbc = (rp_bins[1:]+rp_bins[:-1])/2

In [14]:
print rbc


[  0.12239095   0.17719995   0.25655347   0.37144299   0.53778224
   0.77861137   1.12728837   1.63210958   2.36299934   3.42119545
   4.95327193   7.17144145  10.38294952  15.03263208]

In [15]:
from scipy.stats import normaltest, skewtest, kurtosistest
for i in xrange(wp_vals.shape[1]):
    plt.hist(np.log10(wp_vals[:,i]), bins = int(np.ceil(np.sqrt(1000)) ));
    stat, p = normaltest(np.log10(wp_vals[:,i]) )
    z, ps = skewtest(np.log10(wp_vals[:,i]))
    z, pk = kurtosistest(np.log10(wp_vals[:,i]))
    plt.title(r"$r = %.2f\; $Mpc$ \;\; p=%.2f,%.2f,%.2f $"%(rbc[i], p, ps, pk) )
    plt.show()



In [16]:
sat, p = normaltest(np.log10(wp_vals))
print p


[ 0.87738541  0.32064472  0.77986973  0.68017301  0.6234466   0.1701798
  0.90555643  0.94249975  0.81604299  0.38857586  0.05904191  0.34233148
  0.80526604  0.10194423]

In [ ]: